version 8.0
log using logfilename.log, replace
#delimit ;

*     ***************************************************************** *;
*     ***************************************************************** *;
*       File-Name:      twomodifyingvariables.do                        *;
*       Date:           Month/Day/Year                                  *;
*       Author:         MRG                                             *;
*       Purpose:        Do-file for creating figures showing the        *;
*                       marginal effect of X on Y with two modifying    *;
*                       variables where the dependent variable is       *;
*                       continuous.  The estimated model is:            *;
*                       Y = b0 + b1X  + b2Z  + b3W                      *;
*                              + b4XZ + b5XW + b6ZW + b7XZW + epsilon   *;
*       Input File:     name_of_datafile.dta                            *;
*       Output File:    logfilename.log                                 *;
*       Data Output:    None                                            *;
*       Previous file:  None                                            *;
*       Machine:        Office                                          *;
*     ****************************************************************  *;
*     ****************************************************************  *;


use h:name_of_datafile.dta;

*     ****************************************************************  *;
*     Estimate Model: Y = b0 + b1X  + b2Z  + b3W                        *;
*                            + b4XZ + b5XW + b6ZW + b7XZW + b8Controls  *;
*                            + epsilon                                  *;
*     ****************************************************************  *;

regress Y   X   Z   W   XZ  XW  ZW  XZW Controls;

*     ****************************************************************  *;
*       Generate the values of Z for which you want to calculate the    *;
*       marginal effect (and standard errors) of X on Y.  Z will be     *;
*       the X-axis.                                                     *;
*     ****************************************************************  *;

generate MV=((_n-1)/10);

replace  MV=. if _n>60;

*     ****************************************************************  *;
*       Generate the values of W for which you want to calculate the    *;
*       marginal effect (and standard errors) of X on Y.                *;
*     ****************************************************************  *;

scalar W0=0;

scalar W1=1;

scalar W2=2;

scalar W3=3;

*     ****************************************************************  *;
*       Grab elements of the coefficient and variance-covariance matrix *;
*       that are required to calculate the marginal effect and standard *;
*       errors.                                                         *;
*     ****************************************************************  *;

matrix b=e(b);
matrix V=e(V);

scalar b1=b[1,1];
scalar b2=b[1,2];
scalar b3=b[1,3];
scalar b4=b[1,4];
scalar b5=b[1,5];
scalar b6=b[1,6];
scalar b7=b[1,7];

scalar varb1=V[1,1];
scalar varb2=V[2,2];
scalar varb3=V[3,3];
scalar varb4=V[4,4];
scalar varb5=V[5,5];
scalar varb6=V[6,6];
scalar varb7=V[7,7];

scalar covb1b4=V[1,4];
scalar covb1b5=V[1,5];
scalar covb1b7=V[1,7];
scalar covb4b5=V[4,5];
scalar covb4b7=V[4,7];
scalar covb5b7=V[5,7];

scalar list b1 b2 b3 b4 b5 b6 b7 varb1 varb2 varb3 varb4 varb5 varb6 varb7 
            covb1b4 covb1b5 covb1b7 covb4b5 covb4b7 covb5b7;
            
*     ****************************************************************  *;
*       We want to calculate the marginal effect of X on Y for all      *;
*       MV values of the modifying variable Z. We also want to          *;
*       calculate this marginal effect as Z changes for specific values *;
*       of the second modifying variable W.  In the code below, we      *;
*       calculate the marginal effect of X on Y for all values of Z     *;
*       when W=0, when W=1, when W=2, and when W=3.                     *;
*     ****************************************************************  *;


gen conb0=b1+b4*MV+b5*W0+b7*(MV*W0) if _n<60;

gen conb1=b1+b4*MV+b5*W1+b7*(MV*W1) if _n<60;

gen conb2=b1+b4*MV+b5*W2+b7*(MV*W2) if _n<60;

gen conb3=b1+b4*MV+b5*W3+b7*(MV*W3) if _n<60;


*     ****************************************************************  *;
*       Calculate the standard errors for the marginal effect of X on Y *;
*       for all MV values of the modifying variable Z. Do this for the  *;
*       case when W=0, when W=1, when W=2, and when W=3.                *;
*     ****************************************************************  *;

gen conse0=sqrt(varb1
                + varb4*(MV^2) + varb5*(W0^2) + varb7*(MV^2)*(W0^2)
                + 2*MV*covb1b4 + 2*W0*covb1b5 + 2*MV*W0*covb1b7 + 2*MV*W0*covb4b5
                + 2*W0*(MV^2)*covb4b7 + 2*(W0^2)*MV*covb5b7 if _n<60;

gen conse1=sqrt(varb1
                + varb4*(MV^2) + varb5*(W1^2) + varb7*(MV^2)*(W1^2)
                + 2*MV*covb1b4 + 2*W1*covb1b5 + 2*MV*W1*covb1b7 + 2*MV*W1*covb4b5
                + 2*W1*(MV^2)*covb4b7 + 2*(W1^2)*MV*covb5b7 if _n<60;
                
gen conse2=sqrt(varb1
                + varb4*(MV^2) + varb5*(W2^2) + varb7*(MV^2)*(W2^2)
                + 2*MV*covb1b4 + 2*W2*covb1b5 + 2*MV*W2*covb1b7 + 2*MV*W2*covb4b5
                + 2*W2*(MV^2)*covb4b7 + 2*(W2^2)*MV*covb5b7 if _n<60;
                
gen conse3=sqrt(varb1
                + varb4*(MV^2) + varb5*(W3^2) + varb7*(MV^2)*(W3^2)
                + 2*MV*covb1b4 + 2*W3*covb1b5 + 2*MV*W3*covb1b7 + 2*MV*W3*covb4b5
                + 2*W3*(MV^2)*covb4b7 + 2*(W3^2)*MV*covb5b7 if _n<60;      
                
                                
*     ****************************************************************  *;
*                           Create t statistics                         *;
*     ****************************************************************  *;

gen t0=conb0/conse0;

gen t1=conb1/conse1;

gen t2=conb2/conse2;

gen t3=conb3/conse3;


*     ****************************************************************  *;
*       Generate a `shadow' variable that is equal to the marginal      *;
*       marginal effects that have already been calculated.  These will *;
*       be helpful in placing the stars in the figure as you will see.  *;
*     ****************************************************************  *;

gen consb0=conb0;

gen consb1=conb1;

gen consb2=conb2;

gen consb3=conb3;

*     ****************************************************************  *;
*       Make the `shadow' variable be missing if the t score is not     *;
*       larger than the critical level of significance that you want.   *;
*     ****************************************************************  *;

replace consb0 = . if abs(t0)<1.96;

replace consb1 = . if abs(t1)<1.96;

replace consb2 = . if abs(t2)<1.96;

replace consb3 = . if abs(t3)<1.96;

*     ****************************************************************  *;
*       Generate a string variable called str1 that is designated with  *;
*       a star.                                                         *;
*     ****************************************************************  *;

generate str1 txt="*";

*     ****************************************************************  *;
*       Graph the marginal effect of X on Y across the desired range of *;
*       the modifying variable Z. Do this for when W=0, when W=1, when  *;
*       W=2, and for when W=3.                                          *;
*     ****************************************************************  *;


graph twoway   line conb0  MV,  clpattern(solid) clwidth(thin)
        ||  scatter consb0 MV,  mlabel(txt) msymbol(i) mlabsize(vsmall) mlabgap(-1.0) mlabposition(11)
        ||     line conb1  MV,  clpattern(solid) clwidth(thin)
        ||  scatter consb1 MV,  mlabel(txt) msymbol(i) mlabsize(vsmall) mlabgap(-1.0) mlabposition(11)
        ||     line conb2  MV,  clpattern(solid) clwidth(thin)
        ||  scatter consb2 MV,  mlabel(txt) msymbol(i) mlabsize(vsmall) mlabgap(-1.0) mlabposition(11)
        ||     line conb3  MV,  clpattern(solid) clwidth(thin)
        ||  scatter consb3 MV,  mlabel(txt) msymbol(i) mlabsize(vsmall) mlabgap(-1.0) mlabposition(11)        
        ||  ,   
            xlabel(0 1 2 3, labsize(2.5)) 
            ylabel(-0.5 0 0.5 1 1.5, labsize(2.5))
            yscale(noline)
            xscale(noline)
            legend(off)
            yline(0, lcolor(black))  
            title(Marginal Effect of X on Y As Z and W Change, size(4))
            subtitle(" " "Dependent Variable: Y" " " " ", size(3))
            xtitle(Z, size(3))
            ytitle(Marginal Effect of X, size(3))
            xsca(titlegap(2)) 
            ysca(titlegap(4))
            text(0.9 3.3  "W=3", justification(left) size(2.5))
            text(0.49 3.3 "W=2", justification(left) size(2.5))
            text(0.08 3.3 "W=1", justification(left) size(2.5))
            text(0.33 3.3 "W=0", justification(left) size(2.5))
            text(1.5 1 "* indicates significance at the 95% level", justification(left) size(2.5))
            scheme(s2mono) graphregion(fcolor(white))
            graphregion(margin(r=28));
       

*     ****************************************************************  *;
*                 Figure can be saved in a variety of formats.          *;
*     ****************************************************************  *; 

graph export  h:\figure1.eps, replace;

translate @Graph h:\figure1.wmf;

*     ****************************************************************  *;
*                                   THE END                             *;
*     ****************************************************************  *;


 exit;
